C
C     ..................................................................
C
C        SUBROUTINE MTDS
C
C        PURPOSE
C           MULTIPLY A GENERAL MATRIX A ON THE LEFT OR RIGHT BY
C           INVERSE(T),INVERSE(TRANSPOSE(T)) OR INVERSE(TRANSPOSE(T*T))
C           THE TRIANGULAR MATRIX T IS STORED COLUMNWISE IN COMPRESSED
C           FORM, I.E. UPPER TRIANGULAR PART ONLY.
C
C        USAGE
C           CALL MTDS(A,M,N,T,IOP,IER)
C
C        DESCRIPTION OF PARAMETERS
C           A     - GIVEN GENERAL MATRIX WHITH M ROWS AND N COLUMNS.
C           M     - NUMBER OF ROWS OF MATRIX A
C           N     - NUMBER OF COLUMNS OF MATRIX A
C           T     - GIVEN TRIANGULAR MATRIX STORED COLUMNWISE UPPER
C                   TRIANGULAR PART ONLY. ITS NUMBER OF ROWS AND
C                   COLUMNS K IS IMPLIED BY COMPATIBILITY.
C                   K = M IF IOP IS POSITIVE,
C                   K = N IF IOP IS NEGATIVE.
C                   T OCCUPIES K*(K+1)/2 STORAGE POSITIONS.
C           IOP   - INPUT VARIABLE FOR SELECTION OF OPERATION
C                   IOP = 1 - A IS REPLACED BY INVERSE(T)*A
C                   IOP =-1 - A IS REPLACED BY A*INVERSE(T)
C                   IOP = 2 - A IS REPLACED BY INVERSE(TRANSPOSE(T))*A
C                   IOP =-2 - A IS REPLACED BY A*INVERSE(TRANSPOSE(T))
C                   IOP = 3 - A IS REPLACED BY INVERSE(TRANSPOSE(T)*T)*A
C                   IOP =-3 - A IS REPLACED BY A*INVERSE(TRANSPOSE(T)*T)
C           IER   - RESULTING ERROR PARAMETER
C                   IER =-1 MEANS M AND N ARE NOT BOTH POSITIVE
C                                 AND/OR IOP IS ILLEGAL
C                   IER = 0 MEANS OPERATION WAS SUCCESSFUL
C                   IER = 1 MEANS TRIANGULAR MATRIX T IS SINGULAR
C
C        REMARKS
C           SUBROUTINE MTDS MAY BE USED TO CALCULATE THE SOLUTION OF
C           A SYSTEM OF EQUATIONS WITH SYMMETRIC POSITIVE DEFINITE
C           COEFFICIENT MATRIX. THE FIRST STEP TOWARDS THE SOLUTION
C           IS TRIANGULAR FACTORIZATION BY MEANS OF MFSD, THE SECOND
C           STEP IS APPLICATION OF MTDS.
C           SUBROUTINES MFSD AND MTDS MAY BE USED IN ORDER TO CALCULATE
C           THE PRODUCT TRANSPOSE(A)*INVERSE(B)*A WITH GIVEN SYMMETRIC
C           POSITIVE DEFINITE B AND GIVEN A EFFICIENTLY IN THREE STEPS
C           1) TRIANGULAR FACTORIZATION OF B (B=TRANSPOSE(T)*T)
C           2) MULTIPLICATION OF A ON THE LEFT BY INVERSE(TRANSPOSE(T))
C              A IS REPLACED BY C=INVERSE(TRANSPOSE(T))*A
C           3) CALCULATION OF THE RESULT FORMING TRANSPOSE(C)*C
C
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
C           NONE
C
C        METHOD
C           CALCULATION OF X = INVERSE(T)*A IS DONE USING BACKWARD
C           SUBSTITUTION TO OBTAIN X FROM T*X = A.
C           CALCULATION OF Y = INVERSE(TRANSPOSE(T))*A IS DONE USING
C           FORWARD SUBSTITUTION TO OBTAIN Y FROM TRANSPOSE(T)*Y = A.
C           CALCULATION OF Z = INVERSE(TRANSPOSE(T)*T)*A IS DONE
C           SOLVING FIRST TRANSPOSE(T)*Y = A AND THEN T*Z = Y, IE.
C           USING THE ABOVE TWO STEPS IN REVERSE ORDER
C
C     ..................................................................
C
      SUBROUTINE MTDS(A,M,N,T,IOP,IER)
C
C
      DIMENSION A(1),T(1)
      DOUBLE PRECISION DSUM
C
C        TEST OF DIMENSION
      IF(M)2,2,1
    1 IF(N)2,2,4
C
C        ERROR RETURN IN CASE OF ILLEGAL DIMENSIONS
    2 IER=-1
      RETURN
C
C        ERROR RETURN IN CASE OF SINGULAR MATRIX T
    3 IER=1
      RETURN
C
C        INITIALIZE DIVISION PROCESS
    4 MN=M*N
      MM=M*(M+1)/2
      MM1=M-1
      IER=0
      ICS=M
      IRS=1
      IMEND=M
C
C        TEST SPECIFIED OPERATION
      IF(IOP)5,2,6
    5 MM=N*(N+1)/2
      MM1=N-1
      IRS=M
      ICS=1
      IMEND=MN-M+1
      MN=M
    6 IOPE=MOD(IOP+3,3)
      IF(IABS(IOP)-3)7,7,2
    7 IF(IOPE-1)8,18,8
C
C        INITIALIZE SOLUTION OF TRANSPOSE(T)*X = A
    8 MEND=1
      LLD=IRS
      MSTA=1
      MDEL=1
      MX=1
      LD=1
      LX=0
C
C        TEST FOR NONZERO DIAGONAL TERM IN T
    9 IF(T(MSTA))10,3,10
   10 DO 11 I=MEND,MN,ICS
   11 A(I)=A(I)/DBLE(T(MSTA))
C
C        IS M EQUAL 1
      IF(MM1)2,15,12
   12 DO 14 J=1,MM1
      MSTA=MSTA+MDEL
      MDEL=MDEL+MX
      DO 14 I=MEND,MN,ICS
      DSUM=0.D0
      L=MSTA
      LDX=LD
      LL=I
      DO 13 K=1,J
      DSUM=DSUM-T(L)*A(LL)
      LL=LL+LLD
      L=L+LDX
   13 LDX=LDX+LX
      IF(T(L))14,3,14
   14 A(LL)=(DSUM+A(LL))/T(L)
C
C        TEST END OF OPERATION
   15 IF(IER)16,17,16
   16 IER=0
      RETURN
   17 IF(IOPE)18,18,16
C
C        INITIALIZE SOLUTION OF T*X = A
   18 IER=1
      MEND=IMEND
      MN=M*N
      LLD=-IRS
      MSTA=MM
      MDEL=-1
      MX=0
      LD=-MM1
      LX=1
      GOTO 9
      END
